function trampoline_simulation
    %% === 参数定义 ===
    m = 68;             % 运动员质量 (kg)
    g = 9.81;           % 重力加速度 (m/s^2)
    rho = 1.225;        % 空气密度 (kg/m^3)

    % 飞行阶段参数
    Cd = 1.0;           % 空气阻力系数（姿态相关）
    A = 0.7;            % 迎风面积 (m^2)
    v0 = 5.0;           % 起跳速度 (m/s)
    y0 = 0.93;          % 起跳初始高度 (重心) (m)

    % 蹦床参数
    k = 4951.05;        % 蹦床垂直刚度 (N/m)
    c = 200.0;          % 阻尼系数 (Ns/m)

    %% === 飞行阶段模拟 ===
    fprintf('=== 飞行阶段模拟 ===\n');
    y_init = [y0; v0];             % 初始状态：[y, v]
    tspan = [0 2];                 % 模拟时间范围
    [t1, Y] = ode45(@(t,y) flightODE(t, y, m, g, rho, Cd, A), tspan, y_init);

    % 找落地时刻：y <= 蹦床高度
    trampoline_height = 1.15;
    idx_land = find(Y(:,1) <= trampoline_height, 1);
    if isempty(idx_land)
        error('未检测到落地，请检查初速度或参数');
    end
    v_impact = Y(idx_land, 2);
    fprintf('落地速度 v_impact = %.2f m/s\n', v_impact);

    %% === 落地阶段模拟 ===
    fprintf('\n=== 落地阶段模拟 ===\n');
    x_init = [0; -v_impact];       % 初始状态：[x=0, v=-v_impact]
    [t2, X] = ode45(@(t,x) landingODE(t, x, m, g, k, c), tspan, x_init);

    x_max = max(X(:,1));           % 最大压缩位移
    F_max = k * x_max;             % 最大冲击力
    fprintf('最大压缩位移 = %.4f m\n', x_max);
    fprintf('最大冲击力 = %.2f N\n', F_max);

    %% === 可视化 ===
    figure('Name','蹦床飞行与落地过程','NumberTitle','off');
    subplot(3,1,1);
    plot(t1, Y(:,1), 'b'); ylabel('高度 y(t) [m]');
    title('飞行阶段'); grid on;

    subplot(3,1,2);
    plot(t1, Y(:,2), 'r'); ylabel('速度 v(t) [m/s]');
    title('飞行速度'); grid on;

    subplot(3,1,3);
    plot(t2, k*X(:,1), 'k'); xlabel('时间 t [s]');
    ylabel('冲击力 F(t) [N]'); title('落地冲击力');
    grid on;

end

%% === 飞行阶段微分方程 ===
function dydt = flightODE(~, y, m, g, rho, Cd, A)
    v = y(2);
    Fd = 0.5 * Cd * rho * A * v^2 * sign(v); % 空气阻力
    dydt = [v; -g - Fd/m];
end

%% === 落地阶段微分方程（弹簧-阻尼模型） ===
function dxdt = landingODE(~, x, m, g, k, c)
    dxdt = [x(2); -g - (k/m)*x(1) - (c/m)*x(2)];
end
